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ABSTRACT 


The  well-known,  general-purpose  structural  analysis  program  ICES-STRUDL 
has  been  used  extensively  at  M.I.T.  Lincoln  Laboratory  since  its  public  release 
in  1967.  Since  that  original  release,  several  updates,  issued  by  the  original 
developers  and  the  ICES  Users  Group,  have  produced  changes  in  programs  and 
documentation.  Likewise,  additional  development  and  enhancement  of  ICES-STRUDL 
have  occurred  at  M.I.T.  Lincoln  Laboratory  on  those  areas  (such  as  finite  ele- 
ments, dynamics,  etc.)  specifically  of  interest  to  the  Laboratory. 

This  report  describes  the  finite  element  capabilities  in  the  M.I.T. 

Lincoln  Laboratory  version  of  ICES-STRUDL  and  can  therefore  serve  as  a user’s 
manual.  Sections  are  devoted  to  specification  of  finite  element  geometry,  ele- 
ment properties,  and  element  loading.  Finite  element  modelling  considerations 
are  discussed  and  three  examples  are  presented  to  illustrate  capabilities  often 
needed  at  the  Laboratory,  but  not  described  in  other  ICES-STRUDL  documentation. 
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I.  FINITE  ELEMENT  GEOMETRY  SPECIFICATION 


1,  Planar  Element  Coordinate  System 

Each  planar  finite  element  in  a STRUDL  problem  has  associated  with  it  a 
coordinate  system  in  which  the  results  for  that  element  are  expressed.  This 

coordinate  system  is  independent  of  the  type  of  element  and  depends  only  upon 

the  orientation  of  the  plane  of  the  element.  It  is  called  the  planar  element 
coordinate  system. 

The  use  of  such  a reference  frame  enables  the  element  results  (stresses 
and  strains)  for  all  of  the  elements  in  a particular  plane  to  be  output  with 
respect  to  a common  reference  frame.  This  makes  the  interpretation  of  output 
a relatively  easy  task,  since  similar  results  are  all  expressed  in  the  same 
directions . 

The  planar  element  coordinate  system  of  a triangular  element  (nodes  1,  2, 
3)  is  illustrated  in  Figure  1.1  for  the  general  case. 

The  orientation  of  the  element  affects  the  determination  of  its  planar 

element  coordinate  system.  The  planar  system  is  located  as  follows  (refer  to 
Figure  1.1): 

(1)  The  planar  Zp  axis  is  normal  to  the  plane  of  the  element.  Its 
positive  direction  is  determined  by  applying  the  right-hand  rule 
to  the  node  order  given  in  the  ELEMENT  INCIDENCES  command. 

(2)  The  planar  Xp  axis  coincides  with  the  intersection  of  the  element 

plane  and  the  global  plane. 

* U b 

(3)  The  projection  of  the  planar  axis  onto  the  global  Z axis  is 

r b 

defined  to  be  positive  in  the  same  direction  as  the  global  Z axis. 

b 

(4)  The  planar  Xp  axis  is  defined  to  be  positive  in  the  direction 
determined  by  the  right-hand  rule. 

The  triangular  element  shown  in  Figure  1.1  has  relative  node  numbers  1,  2,  3 
which  correspond  to  the  nodes  listed  in  the  ELEMENT  INCIDENCES  for  the  element. 
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ELEMENT  PLANE 
WITH  XY  PLANE 


Figure  1.1.  Planar  element  coordinate  system  (general). 

The  only  exception  to  the  above  rules  for  location  of  the  planar  element 

coordinate  system  occurs  when  the  element  lies  in  a plane  parallel  to  the 

global  X„Y„  plane.  The  planar  element  coordinate  system  of  such  a triangular 
b b 

element  is  illustrated  in  Figure  1.2. 

The  planar  system  for  an  element  which  lies  in  a plane  parallel  to  the 
global  plane  is  located  as  follows  (refer  to  Figure  1.2): 

(1)  The  planar  Zp  axis  is  normal  to  the  plane  of  the  element.  Its  posi- 
tive direction  is  determined  by  applying  the  right-hand  rule  to  the 
node  order  given  in  the  ELEMENT  INCIDENCES  command. 

(2*)  The  planar  Xp  axis  is  located  parallel  to  the  global  axis  and 

is  positive  in  the  same  direction  as  the  global  X-  axis. 

b 
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(3*)  The  planar  axis  is  defined  to  be  positive  in  the  direction  deter- 
mined by  the  right-hand  rule. 

The  triangular  element  shown  in  Figure  1.2  has  relative  node  numbers  1,  2,  3 
which  correspond  to  the  nodes  listed  in  the  ELEMENT  INCIDENCES  for  the  element. 

2.  LOCAL  and  GLOBAL  Joints 

The  JOINT  COORDINATES  command  has  been  modified  to  provide  the  capability 
for  analyzing  a structure  comprised  of  plates  which  do  not  lie  in  global  planes. 
If  all  finite  elements  attached  to  a joint  are  co-planar  (in  other  than  a 
global  plane),  and  members  are  incident  on  the  joint,  the  joint  could  be 
declared  LOCAL.  This  LOCAL  specification  would  cause  the  solution  at  the  joint 
to  be  carried  out  in  the  planar  element  coordinate  system;  thus  eliminating  po- 
tential instabilities  in  the  stiffness  matrix  caused  by  rotation  into  the 
global  coordinate  system. 
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All  input  for  a LOCAL  joint  must  be  referred  to  the  planar  element  coor- 
dinate system  (e.g.,  JOINT  RELEASES,  JOINT  LOADS,  etc.)*  output  (e.g., 

displacements)  for  a LOCAL  joint  are  also  referred  to  the  planar  element  coor- 
dinate system. 

Joints  connected  only  by  elements  lying  in  one  of  the  three  global  planes 

(X^Y^,  X^Z^,  or  Y^Z^)  should  be  referred  to  as  GLOBAL  joints  for  ease  in  des- 
0 0 0 0 0 0 

cribing  input  and  interpreting  output.  Similarly,  joints  connected  by  co-planar 
elements  (in  other  than  a global  plane)  should  be  referred  to  as  SUPPORTed 
GLOBAL  joints,  and  JOINT  RELEASES  (with  THETA  angles  if  needed)  should  be  used 
to  constrain  the  directions  of  potential  instability. 

Three  cases  exist  in  which  joints  should  not  be  described  as  LOCAL  joints. 
These  are: 

(1)  Nodes  at  which  elements  are  incident  in  more  than  one  plane. 

(2)  Nodes  at  which  incident  elements  are  all  parallel  to  global  planes 
(instabilities  will  not  be  generated  and  the  reduced  number  of 
unknowns  will  automatically  be  used). 

(3)  Nodes  at  which  members  are  incident. 

In  case  (1)  above,  the  specification  of  LOCAL  will  lead  to  incorrect  re- 
sults. In  case  (2),  the  specification  of  LOCAL  will  change  the  orientation 
of  the  input  loads  and  output  results  to  the  planar  element  coordinate  system. 

In  case  (3),  the  node  cannot  be  specified  as  LOCAL,  as  an  error  message  will 
result,  and  it  will  be  processed  as  if  a GLOBAL  specification  had  been  given. 

Joints  7,  8,  9 and  13,  14,  15  in  Figure  1.3  illustrate  case  (1)  above  and 
must  be  described  as  GLOBAL,  since  they  are  nodes  at  which  elements  are  incident 
in  more  than  one  plane.  Nodes  5 and  17  illustrate  case  (2)  and  are  described 
as  GLOBAL  since  their  incident  elements  are  parallel  to  global  planes.  Node  11 
may  be  described  as  a LOCAL  joint  but  is  better  described  as  a SUPPORTed  GLOBAL 
joint  with  associated  JOINT  RELEASES  and  THETA  angles.  Nodes  at  which  the  edge 
members  are  incident  illustrate  case  (3)  and  must  be  described  as  GLOBAL. 
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Edge  Member 


Preferred  Approach; 

JOINT  COORDINATES  $ DEFAULT  IS  GLOBAL 
I TO  10 

11  SUPPORT 

12  TO  21 

UNIT  DEGREES 
JOINT  RELEASES 

11  FORCE  X Y THl  90.  TH3  45. 
Figure  1.3.  Example  of  LOCAL/GLOBAL  joint  designation. 
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3. 


Element  Nodal  Listing 


The  connectivity  of  element  nodes  (or  joints)  is  described  by  listing  the 
joints  associated  with  each  element  in  the  ELEMENT  INCIDENCE  command. 

For  planar  elements  having  only  comer  nodes,  the  joints  are  listed  con- 
secutively around  the  element  boundary,  starting  at  any  arbitrary  corner  node. 
For  planar  elements  having  side  nodes,  the  corner  nodes  are  listed  first  and 
then  the  side  nodes. 

The  mid-side  nodes  for  the  LST  and  LSR  elements  must  be  listed  in  the 
JOINT  COORDINATES  command,  but  one  does  not  have  to  specify  the  joint  coordi- 
nate values.  This  also  pertains  to  the  isoparametric  IPQQ  (or  IPQS)  and  IPCQ 
elements  when  the  boundary  is  straight  and  the  mid-side  node  is  located  at  the 
1/2-  or  1/3-points,  respectively.  The  PRINT  Joint  COORDINATES  command,  given 
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Figure  1.4.  Nodal  sequence  for  planar  elements. 
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after  an  ANALYSIS  command,  will  print  the  computed  joint  coordinates  of  the 
unspecified  mid-side  nodes.  Before  an  ANALYSIS  is  performed,  however,  the 
values  will  be  zero. 

The  nodal  sequence  for  planar  elements  is  illustrated  in  Figure  1.4.  Note 
that  numbering  may  proceed  in  either  a clockwise  or  counterclockwise  direction. 
For  planar  results  to  make  sense  on  a given  plane  however,  it  is  preferable  to 
choose  one  direction  for  all  elements  on  that  given  plane. 

The  relative  node  numbering  for  tridimensional  elements  having  only  corner 
joints,  such  as  the  IPLS  and  TRIP  elements,  is  obtained  by  numbering  an  arbitrary 
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face.  The  face  is  numbered  by  proceeding  around  the  boundary  of  the  face  in  a 
direction  such  that  an  outward-pointing  normal  is  produced  by  applying  the 
right-hand  rule  to  the  node  ordering.  The  back  face  is  then  numbered  by  pro- 
ceeding in  the  same  direction  from  the  joint  directly  back  of  the  first  num- 
bered node.  Refer  to  Figure  1.5. 

The  scheme  for  the  tridimensional  IPOS  element  is  similar.  One  works  from 
the  front  face  to  the  mid-face  nodes  and  then  to  the  back  face.  All  side  nodes 
on  the  front  face  are  numbered  similar  to  planar  element  side  nodes  before  pro- 
ceeding to  the  mid-face  nodes.  Refer  to  Figure  1.6. 


8 


4. 


Convention  for  Stress  Output 


Stress  and  strain  results  for  the  planar  elements  are  output  with  respect 
to  the  element  planar  coordinate  system.  Actual  stresses  are  output  for  the 
plane  stress/strain  elements,  while  stress  resultants  are  output  for  the  plate 
bending  elements.  (These  stress  resultants  are  in  units  of  shear  and  moment 
per  unit  length.)  Refer  to  Figure  1.7  for  positive  sign  convention. 


Figure  1.7".  Sign  convention  for  planar  element  stresses. 

Principal  stress  and  principal  strain  values  are  available  for  some  plane 
stress/strain  elements  with  the  angle  between  the  principal  directions  and  the 
planar  coordinate  system  also  printed. 

Stress  and  strain  results  for  the  tridimensional  elements  are  output  with 
respect  to  the  global  coordinate  system.  Principal  stresses  and  strains  are 
also  available.  Refer  to  Figure  1.8  for  positive  sign  convention. 
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A 


Figure  1.8.  Sign  convention  for  tridimensional  element  stresses. 
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II.  FINITE  ELEMENT  PROPERTY  SPECIFICATION 


1.  Element  Type  Selection 

The  structural  TYPE  command  for  finite  element  analysis  consists  of  options 
for  specifying  various  finite  element  types. 

PLANE  STRESS  and  PLANE  STRAIN  specify  two-dimensional  problems  varying 
only  with  respect  to  Xp  and  Yp.  This  includes  the  surface  and  body  forces, 
which  also  must  be  planar  and  independent  of  Zp,  Two  degrees-of-freedom  are 
defined  at  each  node,  usually  the  u and  v displacements  in  the  element  planar 
coordinate  system. 

In  PLANE  STRESS  problems  the  out-of-plane  normal  stress  is  considered  to 
be  negligible,  while  the  in-plane  stresses  are  assumed  to  be  independent  of  Zp. 
In  other  words, 

a = T = T =0. 
z xz  yz 

and  a , a , T independent  of  z. 

X y xy  ^ 


These  statements  violate  certain  compatibility  equations,  but  can  be  shown  to 
be  a reasonable  assumption  for  thin  plates. 


In  PLANE  STRAIN  problems  the  out-of-plane  strain  is  considered  to  be  zero, 
while  the  in-plane  strains  are  assumed  to  be  independent  of  Zp.  In  other  words. 


and 


e = Y = Y 
z 'xz  'yz 

e , e , Y 

X y 'xy 


0 

independent  of  z. 


The  corresponding  physical  requirement  which  justifies  this  assumption  is  that 
the  thickness  is  significant  relative  to  the  representative  Xp  and  Yp  dimen- 
sions . 

PLATE  BENDING  specifies  problems  where  only  bending  deformations  need  to 
be  considered.  Plate  bending  formulations  include  three  degrees-of-freedom  at 
each  node;  usually  the  displacement  w in  the  direction  normal  to  the  plane 
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of  the  element,  and  the  rotations  0 and  0 in  the  element  planar  coordinate 

X y 

system.  The  effect  of  transverse  shear  flexibility  is  ignored. 

PLATE  specifies  a combination  of  PLANE  STRESS  and  PLATE  BENDING  in  which 
the  uncoupled  stiffnesses  of  the  two  components  are  superimposed  to  represent 
a three-dimensional  thin  flat  plate. 

TRIDIMENSIONAL  specifies  problems  where  a three-dimensional  state  of  stress 
exists  in  the  body  and  requires  solid  elements  rather  than  planar  ones- 

At  least  one  TYPE  command  must  precede  an  ANALYSIS  command-  A suggested 
means  of  specification  is  to  associate  the  TYPE  command  with  the  ELEMENT  INCI- 
DENCES command  and  to  input  the  element  incidence  data  for  all  elements  of  the 
same  type  at  the  same  time.  This  also  makes  it  easy  to  identify  related  ele- 
ments at  a later  time.  Since  the  TYPE  command  is  noted  for  elements  when  pro- 
cessing ELEMENT  INCIDENCES,  it  is  possible  to  combine  PLANE  STRESS  and  PLANE 
STRAIN  elements  in  the  same  ANALYSIS. 

The  currently  available  STRUDL  finite  elements  are  listed  in  Tables  2.1  - 
2.4  for  the  element  types  described  in  this  section.  In  addition,  each  finite 
element  is  described  in  detail  in  the  Appendix.  The  following  properties  are 
available  for  these  elements  unless  noted  otherwise: 

(1)  Isotropic,  orthotropic,  and  anisotropic  material - 

(2)  Lumped  and  consistent  mass  matrices. 

(3)  Element  loadings  (described  in  Section  III)  consisting  of: 

a.  Edge  Forces 

b.  Surface  Forces 

c.  Body  Forces 

d.  Thermal  Forces 

e.  Initial  Strains 
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TABLE  2.1 

PLANE  STRESS/STRAIN  ELEMENTS 


ELEMENT 

NAME 

SHAPE 

NO.  OF 
NODES 

COMMENTS 

CSTG 

3 

PSR 

PSRCSH 

4 t 

1 * 

3 

*2 

4 

IPLQ 

IPLQCSH 

4. 

r^: 

4 

LSI 

1^-4  2 

6 

LSR 

4 

i 

8< 

( 

7 3 

>-  — ■»  i 

> i 

►  • i 

1 5 

»6 

» 

2 

8 

IPQO 

1 5 2 

8 

UTLQl 

4 

IPCO 

4 10  9 ^ 

12 
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TABLE  2.2 

PLATE  BENDING  ELEMENTS 


ELEMENT 

NAME 

SHAPE 

NO.  OF 
NODES 

COMMENTS 

CPT 

3 

No  edge  or  body  forces 

BPR 

4 • • 3 

1 ^ ^2 

4 

No  edge  or  body  forces 

PBQl 

n: 

4 

No  edge  or  body  forces 

TABLE  2.3 
PLATE  ELEMENTS 


ELEMENT 

NAME 

SHAPE 

NO.  OF 
NODES 

COMMENTS 

SBCT 

3 

No  initial  strains 
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TABLE  2.4 

TRIDIMENSIONAL  ELEMENTS 


ELEMENT 

NAME 

SHAPE 

NO.  OF 
NODES 

COMMENTS 

TRIP 

6 

IPLS 

IPLSCSH 

8^ 

4 

I 

1 ^ 

-f7 

1 T 

>2 

8 

IPQS 

5~  * 2 

20 
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2. 


Thermal  Expansion  Coefficients 


The  initial  strains  due  to  temperature  are  evaluated  from: 


At  a 


where 


vector  of  initial  strains 


a 


vector  of  thermal  expansion  coefficients 


AT  = temperature  differential 


The  thermal  expansion  coefficient  directions  correspond  to  the  planar 
element  coordinate  system  for  planar  elements  and  to  the  global  coordinate 
system  for  tridimensional  elements. 

The  ELEMENT  PROPERTIES  command  is  used  to  input  the  thermal  expansion 
coefficients  when  the  material  is  orthotropic  or  anisotropic. 


GTE  [CAX]  v^  [CAY]  v^  [CAZ]  v^ 


[CSXY]  V4  [CSXZ]  v^  [CSYZ]  v^ 

where  ^^2*  ^3  ~ axial  thermal  expansion  coefficients 
V4,  V5,  v^  = shear  thermal  expansion  coefficients 

If  the  material  is  isotropic,  the  shear  thermal  expansion  coefficients 
are  zero  and  the  axial  thermal  expansion  coefficients  are  a constant  value  a. 
We  can  input  a with  the  CONSTANTS  command,  since 
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5>  Rigidity  Matrix  Description 


The  ELEMENT  PROPERTIES  conunand  is  used  to  input  element  type,  thickness, 
and,  optionally,  the  material  rigidity  matrix.  The  planar  stretching  rigidity 
matrix,  is  defined  as  the  relationship  between  in-plane  stresses  and 
strains,  and  is  expressed  as; 


/ a 1 

f £ 

1 ^ 
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The  planar  bending  rigidity  matrix,  relates  the  stress  couples  and  the 
bending  deformations  as  follows: 
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The  rigidity  matrix,  for  the  tridimensional  case  is  defined  in  the  same 

manner  as  the  planar  stretching  rigidity  matrix,  i.e.. 
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Since  the  thickness,  h,  of  the  plate  bending  elements  is  assigned  constant. 


D 


B 


D 


S 


If  the  material  is  orthotropic  or  anisotropic,  we  input  or  using 

the  ELEMENT  PROPERTIES  command.  The  tridimensional  rigidity  matrix  is  input 

by  rows  (6x6  array) , as  is  the  planar  rigidity  matrix  (3x3  array) . The 

-2  ^ 

units  are  expressed  in  FL  terms  (usually  psi) , but  depend  on  the  UNITS  speci- 
fications. 

RIGIDITY  MATRIX  NC  ic 


ROW 


V 


11 


V 


In 


ROW  i V , - V 

n nl  nn 


where : 

ic  = number  of  colvimns  in  the  -rigidity  matrix; 

ij^  . . . i^  = row  numbers  of  the  rows  being  specified; 

V,,  . . V = values  in  the  ridigity  matrix. 

The  rigidity  matrix  for  planar  elements  is  specified  in  the  planar  ele- 
ment coordinate  system,  while  the  rigidity  matrix  for  tridimensional  elements 
is  specified  in  the  global  coordinate  system. 


If  the  material  is  isotropic,  we  input  E,  v with  the  CONSTANTS  command 
since  JD'g  for  plane  strain  can  be  computed  internally  as: 


E 

(1+v)  (l-2v) 


(1-v)  V 0 

V (1-v)  0 

0 .0 
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and  D„  for  plane  stress  as: 
~S  ^ 


D 


■S 


1 V 0 

V 1 0 


0 


0 


1-v 

2 


The  isotropic  tridimensional  rigidity  matrix  may  also  be  expressed  by  in- 
putting E,  V with  the  CONSTANTS  command,  since 


(1+v)  (l-2v) 


(1-v) 

V 

V 
0 
0 
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V 

(1-v) 

V 
0 
0 
0 


V 

V 

(1-v) 

0 

0 

0 
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0 

0 

C1-2V) 

2 

0 


0 

0 

0 


(1-2V) 


0 

0 

0 

0 

0 


(l-2v) 


When  treating  a laminated  plate,  the  equivalent  homogeneous  stretching 


and  bending  ridigity  matrices  are  generated  by  defining  D and  0^  as  follows: 

b D 


i E 

h layers 


D.h. 

~i  1 


E 

layers 


D.dz 


It  sometimes  happens  that  orthotropic  planar  rigidity  exists  in  a local 
element  orientation.  It  is  then  necessary  to  transform  the  orthotropic  local 
rigidity  matrix,  to  an  anisotropic  planar  rigidity  matrix, 


D 


P 
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where  the  stress  transformation^^  is  given  by: 


~a 


cos^6 

• 2. 

sin  9 

sin  29 
2 


sin^9 

2n 

cos  9 

sin  29 
2 


sin  29 
-sin  29 

cos  29 


An  example  of  such  a transformation  is  illustrated  in  Figure  2.1  for  an  ortho- 
tropic edge.  The  positive  sign  convention  for  stresses  and  rotations  in  the 
local  and  planar  reference  frames  is  also  illustrated. 


Figure  2.1.  Example  of  rigidity  matrix  transformation. 
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4,  Stiffness/Mass  Matrix  Description 


The  stiffness  and  mass  matrices  are  expressed  in  partitioned  form,  where 
each  submatrix  is  identified  by  its  hyperrow  and  hypercolumn  subscripts,  i and 
j,  respectively.  Each  submatrix  is  a square  array  whose  order  is  the  number 
of  degrees  of  freedom  per  node.  The  degrees  of  freedom  at  each  node  must  con- 
form to  those  of  the  element  whose  "type  name"  is  given.  The  partitioned  form 
of  the  matrix  is: 
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Since  it  is  symmetric,  only  the  unique,  non-zero  partitions  need  be  given 
(i.e.,  only  the  diagonal  partitions  and  either  the  upper  OR  lower  off-diagonal 
partitions).  If  a submatrix  has  non-zero  terms  only  on  the  diagonal,  the 
DIAGONAL  form  (as  shown  in  the  command  description)  may  be  used. 


MASS 

STIFFNESS 

SUBMATRIX 


(MATRIX)  [NODES]  n [NDF] 


1 D 


v^T  . . . V- 

11  Im 
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V 1 , V 

‘-ml  mm 


DIAGONAL  V,,  ...  V 
\ 11  mm 
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where 


n 


m 

i 


V 


11  ' 


j 

V 

mm 


number  of  nodes  of  the  element; 
number  of  degrees  of  freedom  per  node; 
hyperrow  subscript  of  the  submatrix  being  specified; 
hypercolumn  subscript  of  the  submatrix  being  specified; 
stiffness  values  for  the  submatrix  being  specified. 


The  number  of  hypercolumns  (or  hyperrows)  in  the  element  stiffness/mass 
matrix  will  be  equal  to  the  number  of  nodes,  n.  The  number  of  columns  (or 
rows)  within  each  submatrix  is  equal  to  the  number  of  degrees  of  freedom  per 
joint,  m.  For  example,  a triangular  plane  stress  element  having  three  nodes 
will  have  three  hypercolumns  (and  hyperrows),  while  the  number  of  columns  (or 
rows)  within  each  submatrix  will  equal  two  (corresponding  to  the  u,  v nodal 
degrees  of  freedom  at  each  joint). 

Again,  the  stiffness/mass  matrix  is  specified  in  the  planar  element  coor- 
dinate system  for  planar  elements,  and  in  the  global  coordinate  system  for 
tridimensional  elements.  Units  conversion  is  not  performed  for  stiffness  ma- 
trix components.  They  must  be  given  in  metric  units  of  meters , and  newtons . 

5.  Special  DUMMY  Element  Specification 

When  describing  ELEMENT  PROPERTIES,  it  is  necessary  to  specify  an  element 
"type  name"  (e.g.,  CSTG,  IPLQ,  etc.)  even  if  the  stiffness/mass  matrix  is  in- 
put directly.  This  ^'typre  name"  is  used  to  identify  the  program  name  to  be 
called  within  certain  element  procedures  (e.g.,  stiffness/mass  matrix  genera- 
tion, loads  generation,  stress  backsubstitution,  stress/strain  output,  etc.). 

For  this  reason  a special  type  name,  DUMMY,  has  been  added  to  allow  ele- 
ment information  such  as  stiffness/mass  matrices  to  be  input  without  specify- 
ing a unique  element  "type  name".  Load  generation,  stress  backsubstitution, 
and  output  are  bypassed  for  this  special  element,  since  it  lacks  the  necessary 
information  to  compute  and  store  them. 
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The  specification  of  degrees  of  freedom  for  this  DUMMY  element  is  given 
in  the  ELEMENT  PROPERTIES  command  as: 


DOF  21  2^  M 


These  nodal  degrees  of  freedom  are  global  degrees  of  freedom,  and  must  be  the 
same  at  all  nodes. 

6.  Element  Similarities  Specification 

When  an  element  is  similar  to  another,  it  is  possible  to  take  advantage  of 
element  similarities  to  reduce  the  computation  time.  For  this  purpose,  the 
two  elements  need  not  be  identical.  Differences  may  exist  in  applied  loads, 
boundary  node  support  conditions,  and  geometric  orientation.  This  similarity 
is  specified  with  the  ELEMENT  PROPERTIES  command  as: 


SAME  (AS) 


' PLANAR  [TH4] 

[THI]  [TH2]  [TH3] 


CBUT  LOADINGS] 


where : 


element  name 

values  of  the  global  rotation  angles  0^,  O^,  9^ 
value  of  the  planar  rotation  angle  6^ 


When  an  element  is  identical  to  another,  only  the  identifier  for  the 
second  element  needs  to  be  given  in  the  similarity  specification.  If,  however, 
their  geometric  orientation  with  respect  to  the  global  (or  planar]  coordinate 
system  is  different,  it  is  necessary  to  specify  rotation  angles.  The  angles 
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THl , TH2,  and  TH3  indicate  how  the  base  element  is  to  be  rotated  in  the 
global  reference  frame  so  that  it  coincides  with  the  similar  element  (refer 
to  Figure  2.2) . 


Xg 


Figure  2.2.  Rotation  angles  for  element  similarities. 


Note  that  if  planar  elements  are  to  be  rotated,  the  three-axis  rotation 
procedure  is  extremely  cumbersome.  Since  the  element  planar  coordinate  sys- 
tems are  known  for  both  elements,  it  is  possible  to  use  the  angle  TH4  to  in- 
dicate how  the  planar  base  element  is  to  be  rotated  in  the  planar  reference 
frame  so  that  it  coincides  with  the  similar  planar  element. 
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When  the  element  loads  applied  to  the  similar  element  differ  from  those 
for  the  base  element,  it  is  necessary  to  specify  the  option  BUT  LOADINGS. 
However,  this  is  not  necessary  if  the  only  difference  in  loadings  involves 
applied  joint  loads  or  displacements. 

7.  Examples  of  ELEMENT  PROPERTIES  Specification 

(1)  Element  type  and  thickness 

ELEMENT  PROPERTIES 

1 TO  15  TYPE  'LSR'  THICK  0.5 


(2)  DUMMY  Element  Type  and  Stiffness/Mass  Matrix 

ELEMENT  PROPERTIES 

1 TO  10  TYPE  'DUMMY'  DOF  DX  DZ  - 
STIFFNESS  MATRIX  NODES  2 NDF  2 

SUBMATRIX  1 1 DIAG  87575.  98700. 

SUB  1 2 51500.  47000. 

35000.  63735. 

SUB  2 2 DIAG  75000.  841000. 

MASS  MATRIX  NODES  2 NDF  2 

SUB  1 1 DIAG  .015  .015 

SUB  2 2 DIAG  .015  .015 


(3)  Element  type,  thickness,  thermal  expansion  coefficients,  and 
rigidity  matrix. 

ELEMENT  PROPERTIES 

1001  TO  1016  TYPE  'CSTG'  THICK  0.1  - 

CTE  10.33E-6  10.67E-6  CSXY  0.94E-6 

RIGIDITY  MATRIX  NC  3 

ROW  1 0.15156D+07  0,43313D+06  0.46849D+05 

ROW  2 0.43313D+06  0. 15891D+07  0.57063D+05 

ROW  3 0.46849D+05  0.57063D+05  0.57226D+06 
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(4)  Element  type  and  similarities 
ELEMENT  PROPERTIES 


1 

TYPE 

'PSR' 

' THICK  1. 

2 

TO 

60 

SAME 

AS  1 

61 

TO 

120 

SAME 

AS  1 

ROTATED  TH4  90 
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III.  FINITE  ELEMENT  LOADING  SPECIFICATION 


1.  Edge  Forces 

Element  loads  acting  on  the  exterior  boundaries  of  plane  stress/strain 
finite  elements  or  on  the  edges  of  tridimensional  finite  elements  must  be  con- 
verted into  consistent  joint  loads.  These  consistent  joint  loads  are  defined 
as  *'work  equivalent**  forces: 


nodes 

S P.  • 6u. 
1 1 


6u  dS 


where  u^,  P^  are  the  nodal  displacements  and  equivalent  joint  loads,  respect- 
ively. As  can  be  seen  by  the  right-hand  side  of  the  above  equation,  these 
equivalent  joint  loads  depend  on  the  assumed  element  displacement  expansion, 
u,  as  well  as  on  the  load  intensities,  p. 

EDGE  FORCES  may  be  specified  under  the  ELEMENT  LOADS  command  as: 


list  EDGE  FORCES  EDGE  i^ 


/ * 

LOCAL 

► PLANAR 

LY 

V- 

2 

GLOBAL 

LZ 

v^ 

1 ■ 

^ 1 

\ (variable)  } 

(variable)  = VARIABLE 


< 

X 

Xt 

x„  . 

. . X 

— 

1 

2 

n 

VY 

>^1 

YZ  . 

• • y 

^n 

VZ 

2- 

. 

. . z 

1 

2 

n 

where: 


list 


list  of  elements  for  which  edge  forces  are  applied; 

edge  number  where  load  is  being  applied; 

value  of  the  load  components  for  uniform  loading; 

value  of  X,  Y,  Z load  components  at  node  j for 
variable  loading. 


(Units  of  applied  load  are  force/unit  length) 
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The  edge  numbers  for  planar  elements  are  defined  by  the  nodal  order  in 
the  ELEMENT  INCIDENCES  command.  The  first  pair  of  comer  nodes  define  EDGE  1, 
the  second  pair  EDGE  2,  etc.,  until  all  edges  have  been  defined. 


TN-77-17(3-l) 


EDGE  3 


EDGE  2 


Figure  3.1.  Edge  numbering  for  planar  elements. 

If  the  VARIABLE  loading  specification  is  used,  the  nodal  ordering  on  the  edge 
proceeds  from  comer  node  to  comer  node  in  the  same  direction  as  the  ELEMENT 
INCIDENCE  ordering  (e.g.,  nodes  3,  9,  10,  4 define  edge  3 in  Figure  3.1). 

The  edge  numbers  for  tridimensional  elements  are  also  defined  by  the  ele- 
ment nodal  listing  in  the  ELEMENT  INCIDENCES  command.  The  first  numbered  face 
(front  face)  defines  the  first  group  of  edge  numbers,  the  edges  connecting  the 
front  face  and  back  face  define  the  second  group  of  edge  numbers,  and  the  back 
face  edges  are  defined  last.  Edge  numbering  for  the  IPLS,  IPQS,  and  TRIP  tri- 
dimensional elements  are  illustrated  in  Figures  3. 2-3.4,  respectively. 
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TN-77-17(3-2) 


Figure  3.2.  Edge  numbering  for  IPLS  element. 
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Figure*  5.3.  Edge  numbering  for  IPQS  element. 
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TN-77-17(3-4) 


Edge  Node  ordering 
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Figure  3.4.  Edge  numbering  for  TRIP  element. 


EDGE  FORCES  acting  on  tridimensional  elements  must  be  specified  as 
GLOBAL  loads.  If  EDGE  FORCES  acting  on  planar  elements  are  specified  as 
GLOBAL,  the  loads  will  be  transformed  into  the  PLANAR  element  coordinate  sys- 
tem for  consistent  load  calculations. 

For  planar  stress/strain  elements  a LOCAL  coordinate  system  is  established 
at  each  node  in  the  following  manner  (refer  to  Fig.  3.5): 

(a)  Z.  is  determined  by  nodal  ordering  of  the  ELEMENT  INCIDENCES  command. 

Xj 

(b)  is  the  normal  to  the  edge  at  a joint  (positive  in  the  outward- 
pointing direction) . 

(c)  Y,  is  the  tangent  to  the  edge  at  a joint  (positive  direction  being 

Li 

determined  from  right-hand  rule  on  Z and  X . 

Li  Li 


Figure  3.5.  Local  coordinate  system  for  edge  loads. 
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2. 


Surface  Forces 


Element  loads  acting  on  the  surfaces  of  planar  or  tridimensional  elements 
must  also  be  converted  to  consistent  joint  loads  (refer  to  Edge  Forces  for 
description  of  this  conversion) . 

SURFACE  FORCES  are  specified  under  the  ELEMENT  LOADS  command  as: 


list  SURFACE  FORCES  (FACE  i^) 


LOCAL 

PLANAR 

GLOBAL 


PX  v^ 
PY  v^ 
PZ  V3 


(variable) 


where:  list 


(variable) 


list  of  elements  for  which  surface  forces  are  applied; 

face  number  of  tridimensional  element  where  load 
is  applied; 

value  of  the  load  components  for  uniform  loading; 
refer  to  Edge  Forces  description. 


(Units  of  applied  load  are  force/unit  area) . 


SURFACE  FORCES  can  be  specified  as  either  PLANAR  or  GLOBAL  for  planar 
finite  elements.  The  LOCAL  option  is  not  allowed.  If  SURFACE  FORCES  acting 
on  planar  elements  are  specified  as  GLOBAL,  the  loads  will  be  transformed  into 
the  planar  element  coordinate  system  for  consistent  load  calculations.  The 
FACE  i specification  must  not  be  given  for  planar  elements,  since  it  is  clear 
where  the  load  is  to  be  applied.  The  variable  option  is  not  available  for 
plate  bending  elements. 

SURFACE  FORCES  can  be  specified  as  either  LOCAL  or  GLOBAL  for  tridimen- 
sional elements.  The  LOCAL  coordinate  system  refers  to  the  transformed  curvi- 
linear coordinates  (refer  to  Figure  3.6). 

(a)  Z is  the  outward-pointing  normal  to  the  surface, 

Lj 

(b)  X is  parallel  to  edges  1 and  3 of  the  transformed  coordinate  axes. 

Li 

(c)  is  parallel  to  edges  2 and  4 of  the  transformed  coordinate  axes. 


31 


CURVILINEAR  COORDINATES  TRANSFORMED  LOCAL 

COORDINATES 


Figure  3.6.  Local  coordinate  systems  for  surface  loads. 

The  face  numbers  for  tridimensional  elements  are  defined  by  the  nodal 
order  in  the  ELEMENT  INCIDENCES  command.  The  first  numbered  face  is  always 
FACE  1,  and  the  back  face  is  always  FACE  2.  The  connecting  faces  are  then 

numbered  in  the  same  direction  as  the  node  numbering  (FACE  3 ►FACE  N)  . 

The  nodal  ordering  for  each  face  is  started  at  the  first  node  encountered  for 
that  face  in  the  nodal  list.  The  remaining  nodes  on  the  face  are  ordered  by 
proceeding  around  the  face  in  a direction  such  that  an  out vard -pointing  normal 
is  produced  by  applying  the  right-hand  rule  to  the  node  ordering. 
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Figure  3.7.  Face  numbering  for  TRIP  element. 
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Figure  3.8.  Face  numbering  for  IPLS  element. 
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Figure  3.9.  Face  numbering  for  IPQS  element. 


3.  Body  Forces 

Forces,  such  as  dead  load,  often  act  on  a finite  element  throughout  the 
entire  body  of  the  element.  These  body  forces  must  be  converted  into  consis- 
tent joint  loads  (as  detailed  under  Edge  Forces)  for  each  finite  element. 

BODY  FORCES  may  be  specified  under  the  ELEMENT  LOADS  command  as: 


where : 


list 


f * 

^1 1 

i — ►planar)  I 

BODY  FORCES 

) GLOBAL ( i 

^2  1 

BZ 

1 / 

^7  1 

3 I 

(variable) 

list 


Vi,  v^,  Vj 


(variable) 


list  of  elements  for  which  body  forces  are  prescribed; 
value  of  the  load  components  for  uniform  loading; 

refer  to  Edge  Forces  description. 


(Units  of  applied  load  are  force/unit  volume) . 
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BODY  FORCES  must  be  specified  only  as  GLOBAL  for  tridimensional  elements. 

However,  they  may  be  specified  as  either  PLANAR  or  GLOBAL  for  planar  finite 
elements.  If  BODY  FORCES  for  planar  finite  elements  are  specified  as  GLOBAL, 
the  loads  will  be  transformed  into  the  planar  element  coordinate  system  for 
consistent  load  calculation. 

4.  Thermal  Forces 

Temperature  changes  and  thermal  gradients  acting  on  finite  elements  also 
require  conversion  to  equivalent  consistent  joint  loads  (as  detailed  under 
Edge  Forces)  for  each  finite  element. 

Thermal  forces  may  be  specified  under  either  the  ELEMENT  or  JOINT  TEMPERA- 
TURE  commands  as : 


f CHANGE 

1 * 

1 

X 

^2 

1 — ► ipCAL 

list  i 

1 GRADIENT 

Y 

^3 

^^4 

/ GLOBAL 

where : 

list 

= a list  of  joints/elements  at  which 
forces  are  applied; 

the  thermal 

^1 

= value  of  temperature  change,  T; 

V2,  V3, 

= value  of  X,  Y,  Z thermal  gradients 
9T/3Y,  and  9T/3Z) . 

(i.e.,  3T/3X 

Coefficients  of  thermal  expansion  must  be  given  for  all  elements  connected  to 
joints  which  are  subjected  to  JOINT  TEMPERATURES  and  for  all  elements  subjected 
to  ELEMENT  TEMPERATURES. 

When  only  a coarse  description  of  the  thermal  distribution  is  known,  it  is 
often  simplest  to  use  the  ELEMENT  TEMPERATURE  form  of  the  command.  On  the 
other  hand,  if  variations  in  temperature  are  well  defined  and  rapidly  changing 
across  an  element,  it  would  be  better  to  use  the  JOINT  TEMPERATURE  form  of  the 
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command.  In  this  latter  case,  equivalent  nodal  forces  due  to  temperature  are 
generated  by  assuming  the  temperature  expansion  coincides  with  the  displacement 
expansion. 

Plane  stress/strain  elements  and  tridimensional  elements  use  only  the 
temperature  CHANGE  to  generate  thermal  forces,  while  plate  bending  elements 
use  only  temperature  GRADIENT  information.  Only  the  plate  elements  make  use 
of  both  temperature  CHANGE  and  temperature  GRADIENT  data. 

Temperature  GRADIENTS  can  be  either  LOCAL  or  GLOBAL.  However,  if  LOCAL 
is  specified,  the  only  applicable  load  component  is  Z which  defines  a gradient 
normal  to  the  surface.  If  GLOBAL  is  given,  the  component  of  the  gradient  vec- 
tor normal  to  the  plate  will  be  calculated  and  used.  The  in-place  components 
will  be  discarded. 


5.  Initial  Strains 

If  any  initial  strains  are  present  in  a finite  element  from  such  sources 
as  imperfect  fit,  etc.,  then  they  must  be  accounted  for  in  the  model.  The 
elastic  strain  (i.e.,  the  only  one  that  produces  stress)  can  be  expressed  as 
e = e - - Gj,  where  e is  the  total  strain,  is  the  thermal  strain,  and 

Gj  is  the  initial  strain. 

Since  the  deformed  shape  of  the  structure  under  total  loading  (external 
loads,  thermal  strain,  initial  strain)  is  desired,  equivalent  consistent  joint 
loads  are  generated  for  the  thermal  strain  and  initial  strain  effects.  The 
final  elastic  strains  (those  computed  by  STRUDL)  are  then  obtained  by  sub- 
tracting the  input  values  of  thermal  strain  and  initial  strain  from  the  total 
strains  computed  from  the  displacements. 

These  initial  strains  may  be  specified  with  the  ELEMENT  INITIAL  STRAINS 
command  as : 


list 


AXIAL  [X]  [Y^]  V2  [Z] 

SHEAR  [2^]  [2^]  [;^] 
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where : 


list 


a list  of  elements  for  which  initial  strains  are 
prescribed; 

values  of  initial  axial  strains; 


values  of  initial  shear  strains. 


The  initial  strain  components  must  be  given  in  the  GLOBAL  reference  frame 
for  tridimensional  elements  and  in  the  PLANAR  reference  frame  for  planar 
elements . The  strain  curvatures,  rather  than  strains,  are  specified  in  the 
case  of  plate  bending  elements.  Therefore,  initial  strains  must  not  be 
specified  for  plate  elements,  as  inconsistencies  will  result.  (Plate  elements 
may  be  remodeled  as  overlaid  plane  stress  and  plate  bending  elements  if  ini- 
tial strains  are  present.) 
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IV.  FINITE  ELEMENT  MODELING  CONSIDERATIONS 


Convergence  and  Accuracy  of  Model 

One  objective  of  the  analyst  when  using  finite  elements  to  model  a struc- 
ture is  to  generate  results  sufficiently  accurate  to  satisfy  his  requirements. 
To  have  confidence  in  his  results,  he  must  know  their  level  of  accuracy  (i.e., 
establish  a ’^percentage  of  error”  from  the  ”true  answers”) . He  also  wants  his 
results  to  converge  to  the  ”true  solution”  as  the  finite  element  mesh  is  made 
finer. 

Two  basic  criteria  must  be  satisfied  to  assure  convergence. 

(1)  Any  required  state  of  constant  strain  (including  rigid  body 
movement)  can  be  reproduced  on  an  element; 

(2)  Displacements  must  be  continuous  between  adjacent  elements. 

If  these  criteria  are  met,  the  solution  represents  an  upper -bound  on  the  total 
potential  energy.  The  potential  energy  can  be  made  to  converge  monotonically 
if  successively  finer  modeling  retains  the  nodes  of  the  previous  models.  If 
the  second  criterian  is  unsatisfied,  these  non-compatible  elements  will  still 
converge  if  the  first  criterion  above  is  satisfied.  However,  a monotonic  con- 
vergence cannot  be  guaranteed,  and  bounds  on  the  potential  energy  are  not  es- 
tablished. Non-compatibility  is  often  associated  with  plate  bending  elements 
(e.g.,  CPT  and  SECT).  However,  these  non-compatible  elements  have  often  been 
shown  to  exhibit  better  convergence  characteristics  than  compatible  elements 
which  use  complex  displacement  expansions. 

In  generating  the  finer  mesh,  nodes  may  be  added,  but  none  may  be  removed 
or  eliminated.  In  addition,  the  displacement  expansions  of  the  elements  must 
not  be  a function  of  the  element  size  or  orientation. 

Since  convergence  can  be  assured  for  all  finite  element  models,  one  sim- 
ple method  of  assessing  accuracy  is  available.  The  problem  can  be  analyzed 
with  successively  finer  meshes,  whose  successive  results  will  approach  the 
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true  solution.  The  results  of  successive  refinements  as  a function  of  the 
number  of  degrees  of  freedom  allow  the  true  solution  to  be  approximated  by 
extrapolation.  This  process  of  successive  refinements  can  be  terminated 
whenever  the  analyst  feels  he  has  achieved  sufficient  accuracy. 

It  must  be  remembered  that  the  degree  of  refinement  required  for  a parti> 
cular  analysis  is  a function  of  both  the  geometry  of  the  structure  being  ana- 
lyzed and  its  loading  environment.  Finer  meshes,  or  higher-order  elements, 
are  generally  required  in  areas  of  high  stress  gradients.  High  stress  gradi- 
ents are  indicated  by  large  stress  differences  in  adjacent  elements.  If 
structural  loadings  do  not  generate  any  such  high  stress  gradients,  a finer 
grid  will  not  be  required.  Similarly,  if  a curved  boundary  is  being  modeled, 
elements  with  straight  edges  can  never  model  the  problem  exactly,  but  several 
higher-order  elements  may  give  excellent  results  with  significantly  less  ex- 
pense . 

Convergence  and  accuracy  studies  should  be  standard  procedure  for  any 
finite  element  analysis  in  order  that  the  analyst  can  estimate  his  errors  and 
account  for  them. 

2.  Interpretation  of  Results 

To  obtain  the  results  of  a finite  element  analysis,  it  is  necessary  to 
use  the  LIST  command.  Available  for  output  are  the  nodal  displacements  and 
element  stresses  and  strains.  Reactions  and  loads  (statics  check  results) 
are  not  available  in  problems  containing  finite  elements. 

Since  the  stresses  in  a finite  element  problem  are  a function  of  deriva- 
tives of  the  displacements,  the  displacement  results  will  converge  more  quickl 
than  the  stresses.  Therefore,  acceptable  displacement  results  are  often  ob- 
tained for  a mesh  which  is  significantly  coarser  than  that  required  to  obtain 
acceptable  stress  results.  One  technique  which  enables  the  analyst  to  improve 
his  stress  results  without  resorting  to  re-analysis  using  a finer  discretiza- 
tion is  stress  averaging.  The  basic  assumption  is  that  if  the  displacements 
are  accurate,  the  stresses,  on  the  average,  should  also  be  good. 
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Several  variations  of  the  stress  averaging  technique  exist,  all  of  which 
result  in  a smoothing  of  the  stress  pattern  over  the  structure: 

(1)  Nodal  Stress  Averaging 

The  analyst  may  be  interested  in  a set  of  stress  values  at  the  nodes. 
If  this  is  the  case,  the  stresses  from  all  of  the  elements  incident  on 
that  node  may  be  averaged  (the  sum  of  the  stresses  in  a given  direction 
divided  by  the  number  of  incident  elements).  If  the  elements  have  stress 
values  at  the  nodes,  those  stresses  may  be  used.  If  the  elements  have 
constant  stresses,  those  constant  values  may  be  used.  In  either  case, 
the  result  is  a smoothing  of  the  set  of  stress  values  at  particular 
points  in  the  structure. 

(2)  Element  Stress  Averaging 

Alternatively,  the  analyst  may  be  interested  in  a good  approximation 
of  the  stresses  in  the  elements.  If  this  is  the  case,  the  stresses  from 
each  of  the  nodes  may  be  averaged.  If  the  elements  are  constant  stress 
elements  (i.e.,  there  are  no  nodal  stresses),  the  constant  stress  value 
may  be  used. 

It  must  be  noted  that  this  sort  of  averaging  technique  may  not  be  valid 
in  stress  fields  where  the  smoothing  would  tend  to  mask  areas  of  high  stress 
gradients. 

In  summary,  the  intent  of  this  section  has  been  to  point  out  the  fact 
that  while  displacement  results  may  be  sufficiently  accurate  for  a particular 
problem,  it  may  be  appropriate  to  interpret  the  resulting  stresses  somewhat 
more  loosely,  due  to  the  inaccuracies  which  may  have  been  introduced  by  the 
original  finite  element  modeling  assumptions. 

The  output  currently  available  for  the  STRUDL  finite  elements  is  listed 
in  Table  4.1.  In  addition,  each  finite  element  is  also  further  described  in 
the  Appendix. 
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OUTPUT  SUMMARY  FOR  FINITE  ELEMENTS 
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Location  of  output  for  these  elements  can  be  requested  at  edges  or  faces  by  specifying 
DUMP  ON  16384  previous  to  an  ANALYSIS  command. 


V.  EXAMPLES 


1.  Plane  Stress  Example 

A solid  cantilever  beam  (see  Figure  5ol)  will  be  modeled  using  plane  stress 
finite  elements.  The  cross-section  is  rectangular,  12  inches  in  depth,  6 inches 
in  thickness,  and  120  inches  in  length.  A concentrated  end  load  (P  = 10  kips) 
and  a uniform  loading  (w  = 0.8  kips/ft)  are  to  be  investigated.  Tip  deflection 
of  the  beam  is  desired. 


Lj 

P = 10  kips 
W = 0.8  kips /ft  1 

t = 6" 

^ E - 10,000  ksi 

i'  = 0.  0 

* + V + * * ■ 

II 

n 

L = 120" ► 

TN-77-17(5-l) 

Figure  5.1.  Solid  cantilever  beam. 


Since  the  tip  deflection  is  the  desired  output  (rather  than  stresses), 
it  is  possible  to  use  the  CSH  elements  and  reduce  the  number  of  linear  dis- 
placement finite  elements  (PSR,  IPLQ)  required  for  this  model.  It  should  also 
be  noted  that  it  is  possible  to  use  only  one  quadratic  displacement  finite  ele- 
ment (e.g.,  LSR,  IPQQ)  since  the  theoretical  displacement  function  for  this 
example  happens  to  be  quadratic. 

IPLQ  and  IPLQCSH  elements  were  selected  for  the  model.  They  also  were 
used  to  illustrate  how  edge  loads  (w  = 0.8  kips/ft)  can  be  applied  to  the  ele- 
ments. A typical  2-element  model  is  shown  in  Figure  5.2. 


Figure  5.2. 


Finite  element  model  of  beam. 


TN-77-17(5-2) 
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It  should  be  noted  that  adding  joints  7,  8,  9 in  Figure  5.2  does  not  im- 
prove the  tip  displacements  at  all.  Even  though  the  number  of  finite  elements 
is  doubled,  the  element  aspect  ratio  has  also  been  doubled,  and  no  improvement 
in  tip  deflection  is  obtained.  This  illustrates  the  importance  of  proper  dis- 
cretization for  proving  convergence  of  a finite  element  model. 

A typical  input  deck  for  this  example  is  listed  in  Figure  5.3,  while  the 
results  are  illustrated  in  Figure  5.4.  Two  cases  were  considered:  (1)  a long 
beam  (h/L  = 0.1)  in  which  bending  dominates;  and  (2)  a short  beam  (h/L  = 0.5) 
in  which  shear  deformation  is  important. 

It  should  be  noted  that  the  CSH  elements  provide  almost  exact  results  for 
the  loading  cases  illustrated  in  this  example. 

STRUDL  'IPLQCSH'  'PLANE  STRESS  EXAMPLE' 

JOINT  COORDINATES 

1 0.  0.  SUPPORT 

2 0.  12.  SUPPORT 

3 60.  0. 

4 60.  12. 

5 120.  0. 

6 120.  12. 

JOINT  RELEASES 
2 FORCE  Y 

TYPE  PLANE  STRESS 
ELEMENT  INCIDENCES 
113  4 2 
2 3 5 6 4 

ELEMENT  PROPERTIES 
1 2 TYPE  'IPLQCSH'  THICK  6. 

CONSTANTS 
E 10. E6  ALL 
POISSON  0.  ALL 
UNITS  KIPS  FEET 
LOADING  1 
JOINT  LOADS 
6 FORCE  Y -10. 

LOADING  2 
ELEMENT  LOADS 

1 2 EDGE  FORCES  EDGE  3 LOCAL  UNIFORM  LY  -.8 
DUMP  TIME 

STIFFNESS  ANALYSIS 
UNITS  INCHES 

LIST  DISPLACEMENTS  JOINT  5 ItN-77-17 (5-3) 

FINISH  I L_ 

Figure  5.3.  Input  listing  of  STRUDL  model  of  beam. 
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V*  Of  Thaoreticoi  Tip  Deflection 


TN-77-17(5-4) 


CSH 


N N 


Figure  5.4.  Results  of  solid  cantilever  example. 
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2. 


Fracture  Mechanics  Example 


The  quadratic  isoparametric  elements  CIPQQ>  IPQS)  existing  in  STRUDL  have 
been  shown to  embody  singularity  for  calculating  stress  intensity  fac- 

tors of  elastic  fracture  mechanics.  The  singularity  is  obtained  by  placing 
those  mid-side  nodes  near  the  crack  tip  at  the  quarter  points. 

The  STRUDL  input  file  need  not  be  changed  to  accommodate  fracture  analy- 
sis, except  that  the  coordinates  of  the  quarter-side  nodes  must  be  specified. 
The  distance  from  the  quarter-side  nodes  to  the  crack  tip  should  probably  be 
between  2 - 3%  of  the  crack  length.  Elastic  fracture  mechanics  problems  that 
can  be  treated  are  Mode  I (crack  opening),  Mode  II  (crack  sliding),  or  combined 
Modes  I and  II  (refer  to  Figure  5.5). 


TN-77-17(5-5) 


I Opening 


H Sliding 


Figure  5.5.  Elastic  fracture  mechanics  modes. 
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Since  STRUDL  permits  thermal  loading,  the  effects  of  cracks  on  structures  sub- 
jected to  various  thermal  environments  can  also  be  studied. 

In  order  to  compute  the  stress  intensity  factors  and  the  displace- 

^ . \2^ 

ments  at  the  quarter-side  nodes  must  be  found.  The  following  equations  may 
then  be  used  to  evaluate  them. 


u = 


1 ( To  e 

— ^ — /cosa  (2k  - 1)  cos cos  — 

G 1 I'n  \ L 2 2J 

r(2K  + 1)  sin  — - sin  — | + Icosa  r(2K  + 3)  sin  — + sin  — 1 

L 2 2 J ) 4G  T 27T  ( L 2 2J 

To  6 361) 
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|^(2k  + 


6 

(2k  - 1)  cos  — - cos 
2 


7] 


+ cosa 


e 


1)  sin  “ - si 
2 


in  — 1 -4/^  |sina  r(2K  + 3)  si 

2J  ) 4G  f 27T  ( L 


0 . 36 

sin  — + sin  — 


iz-o  T^  6 36 

- cosa  (2k  - 3)  cos  — + cos  — 


where  u^  and  v^  are  the  unknown  displacements  of  the  crack  tip,  u and  v are 
the  displacements  at  the  quarter-side  nodes,  G is  the  shear  modulus,  a is  the 
angle  of  the  crack  as  measured  positive  counterclockwise  from  the  x axis,  and 
r and  6 are  the  crack  tip  centered  polar  coordinates  of  the  quarter-side  nodes. 
The  quantity  k is  dependent  on  the  problem  type: 

f3  - 4V  (plane  strain  or  axisymmetric) 

(3  - V)/(l  + v)  (plane  stress) 
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where  v is  Poisson*s  ratio.  It  should  be  noted  that  this  displacement  assump- 
tion, when  properly  differentiated  to  yield  strains  and  stresses,  yields  the 

-1/2 

correct  singularity  for  stresses  at  the  crack  tip,  i.e.,  (r  ,6) 

which  tends  to  infinity  as  r tends  to  zero. 

Consider  a steel  compact  tensile  specimen,  5 inches  square  with  a 1.5-inch 
crack  and  a 1 psi  tensile  loading  on  the  top  surface.  A 2-element  model  of 
the  specimen  is  shown  in  Figure  5.6. 


1 psi 


TN-77-17(5-6) 


E = 30(10)  psi 
u = 0.  3 
t = 1 inch 


Figure  5.6.  Finite  element  model  of  fracture  specimen. 


The  singularity  of  strains  at  the  crack  tip  is  evidenced  by  the  Jacobian 
of  transformation  from  curvilinear  space  to  Cartesian  space.  This  singularity 
is  bypassed  during  strain  calculations  and  a warning  is  printed  indicating  the 
node  at  which  the  singularity  occurs  (i.e.,  this  should  be  the  node  at  the 
crack  tip).  For  example,  for  the  above  problem 
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****  STRUDL  WARNING  4.35  - JACOBIAN  VANISHED  FOR  ISOPARAMETRIC  ELEMENT  I 
AT  INTEGRATION  POINT  2 

****  STRUDL  WARNING  4.35  - JACOBIAN  VANISHED  FOR  ISOPARAMETRIC  ELEMENT  2 
AT  INTEGRATION  POINT  I 


where  node  3 is  the  integration  point  noted  in  both  warning  messages. 


In  order  to  achieve  enough  accuracy  in  the  output  displacements,  it  is 
necessary  to  apply  a larger  (e.g.,  1000  psi)  tensile  stress  on  the  top  surface 
and  also  request  9 digits  of  output  data  by  the  OUTPUT  DECIMAL  9 command.  The 
input  listing  of  the  data  deck  for  this  example  is  listed  in  Figure  5.7. 


STRUDL  'IPQQ'  'FRACTURE  MECHANICS  EXAMPLE' 
JOINT  COORDINATES 
1 0.  0. 

2 1.125  0. 

3 1.5  0.  SUPPORT 

4 2.375  0.  SUPPORT 
55.  0.  SUPPORT 

6 1.5  .625 

7 

8 

9 0.  2.5 

10 

11  1.5  2.5 

12 

13  5.  2.5 
JOINT  RELEASES 
4 5 FORCE  X 
TYPE  PLANE  STRAIN 
ELEMENT  INCIDENCES 

1 1 3 11  9 2 6 10  7 

2 3 5 13  11  4 8 12  6 
ELEMENT  PROPERTIES 

1 2 TYPE  'IPQQ'  THICK  1. 

CONSTANTS 
E 30. E6  ALL 
POISSON  .3  ALL 
LOADING  1 
ELEMENT  LOADS 

1 2 EDGE  FORCES  EDGE  3 GLOBAL  UNIFORM  LY  - 
DUMP  TIME 

STIFFNESS  ANALYSIS 
OUTPUT  DECIMAL  9 

LIST  DISPLACEMENTS  JOINTS  246 
FINISH 


1000. 


TN-77-17(5-7) 


Figure  5.7.  Input  listing  of  STRUDL  model  of  fracture  specimen. 
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The  number  of  quarter-side  nodes  surrounding  a crack  tip  affects  the  ac- 
curacy of  computation  of  the  stress  intensity  factors.  Several  grids  were  used 

for  the  above  problem  and  convergence  to  the  exact  solution  was  extremely  good. 

1/2 

The  exact  answer  for  Kj  is  4.016  psi-in  . Convergence  is  illustrated  in 
Figure  5.8. 


NUMBER  OF  QUARTER-SIDE  NODES 


Figure  5.8.  Convergence  of  Kj  for  fracture  specimen. 


49 


3.  Super  Finite  Element  Generation  Example 


The  construction  of  super  finite  elements  (or  substructures)  is  often  re- 
quired to  analyze  large  structures.  These  super  finite  elements  contain  only 
those  degrees  of  freedom  needed  for  structural  coupling  and  those  additional 
degrees  of  freedom  needed  to  accurately  reflect  the  stiffness  and  mass  proper- 
ties of  the  structural  component. 

Consider  a flat  plate,  support  frame,  and  torque  tube  assembly  (modeled 
as  simply  as  possible  in  Figure  5.9  for  illustrative  purposes).  The  honeycomb 
plate  section  consists  of  .05"  facesheets  on  a .95"  thick  core. 


50 


The  assembly  connects  to  the  total  structure  at  joints  10  and  12.  Thus,  it  is 
possible  to  obtain  a proper  stiffness  representation  for  the  above  assembly 
with  one  super  finite  element  (or  member)  connecting  only  those  two  joints. 

A proper  mass  representation  would  probably  require  some  additional  degrees 
of  freedom  to  be  considered  for  both  mass  and  stiffness  matrix  development. 

Since  STRUDL  does  not  compute  reactions  at  supports  for  problems  contain- 
ing finite  elements,  it  is  necessary  to  obtain  them  in  another  way.  Two  addi- 
tional joints  *10A’  and  *12A'  are  constructed  differing  only  from  joints  10 
and  12  by  some  very  small  AX  (in  this  case,  .01  inches).  These  two  additional 

o 

joints  are  then  made  SUPPORTS  and  rigid-linked  (accomplished  by  specifying  very 
large  cross-sectional  properties)  by  two  members  to  their  corresponding  real 
joints.  The  MEMBER  FORCES  in  these  two  rigid-link  members  should  be  identical 
to  the  reactions  if  one  is  careful  to  keep  the  sign  convention  straight.  (It 
should  be  simple  to  transform  the  results  from  the  member  coordinate  system 
into  the  global  coordinate  system) . 

The  calculation  of  the  equivalent  plate  thicknesses  to  use  for  the  planar 
ELEMENT  PROPERTIES  is  rather  straightforward  (refer  to  Section  II-3).  The 
modulus  of  elasticity  of  the  plane  stress  and  plate  bending  elements  is  chosen 
to  be  that  of  the  honeycomb  plate.  Their  thicknesses  are  then  computed  to  be: 

(stretching)  h = 2 h.  = 0.10  inch 

(bending)  h^  = 24  h^  ( — - 0.67  inch 

One  interesting  observation  should  be  made  at  this  time.  The  flat  plate 
is  intentionally  modeled  with  triangular  finite  elements  to  illustrate  that 
skew  symmetric  answers  will  be  obtained  for  this  otherwise  perfectly  symmetric 
problem.  This  skew  symmetry  can  be  eliminated  by  overlaying  an  equal  number  of 
finite  elements  whose  diagonals  run  the  other  way,  and  then  modifying  the  plate 
thicknesses  of  all  the  elements  (stretching  h^  = *5  h^  and  bending 
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The  stiffness  matrix  relates  loads  and  displacements.  One  column  of  this 
matrix  can  be  constructed  by  displacing  a SUPPORTed  degree  of  freedom  one  unit, 
while  constraining  all  others,  and  then  listing  the  reactions  at  those  other 
SUPPORTed  degrees  of  freedom.  There  are  12  degrees  of  freedom  needed  for  con- 
struction of  the  super  finite  element  stiffness  matrix,  and  thus,  12  loading 
conditions  are  required. 

Loads  acting  on  the  assembly  are  also  transferred  to  the  coupling  nodes. 
The  reactions  at  the  fully  SUPPORTed  coupling  nodes  represent  the  super  finite 
element  loads  to  be  applied  to  the  coupling  nodes  in  subsequent  models. 

Note  that  OUTPUT  DECIMAL  2 is  specified  in  Figure  5.10  before  listing 

7 

member  forces,  so  that  the  large  stiffness  matrix  values  (of  order  10  ) will 
not  overflow  the  output  field.  Also,  OUTPUT  BY  MEMBER  has  been  specified  so 
that  the  output  from  the  12  loading  cases  will  be  presented  under  each  member 
separately. 

As  mentioned  earlier,  note  the  lack  of  agreement  in  the  marked  terms  in 
Figure  5.11  between  the  skewed  triangular  model  and  the  overlaid  symmetric 
triangular  model.  Non- zero  coupling  terms  are  generated  for  the  skewed  tri- 
angular model  even  though  coupling  does  not  really  exist.  This  illustrates 
very  vividly  that  care  should  be  exercised  when  constructing  finite  element 
models  with  triangular  elements. 

The  member  forces  (really  reactions)  are  then  used  as  stiffness  matrix 
data  for  input  to  another  model  with  n^  loss  of  information  by  this  reduction 
in  problem  size.  While  a MEMBER  could  have  been  selected  to  input  the  stiff- 
ness matrix  data  just  generated,  a DUMMY  finite  element  was  chosen  instead, 
in  order  to  illustrate  the  input  technique  for  finite  element  stiffness  matrix 
input.  The  input  portion  of  the  data  deck  pertaining  to  the  DUMMY  element 
specification  is  listed  in  Figure  5.12. 


52 


STBODL  'SOPEB*  •SDPEH  FIMITE 
JOIMT  COOITDIHATES 

1 0.  0.  15. 

2 25.  0.  15. 

3 50.  0.  15. 

4 0.  0.  0. 

5 25.  0.  0. 

6 50-  0.  0. 

7 0.  0.  -15. 

8 25-  0.  -15. 

9 50.  0.  -15. 

10  0.  -10.  2. 

MOA*  -.01  -10.  2.  SOPPOBT 

11  25.  -10.  2. 

12  50.  -10.  2. 

M2A'  50.01  -10.  2.  SOPPOBT 
TYPE  PLAME  STRESS 
ELEBBNT  IfCIDEKCES 

101  421 

102  245 

103  152 

104  5 1 4 

105  5 3 2 

106  356 

107  263 

108  625 

109  754 

110  578 

111  485 

112  847 

113  865 

114  689 

115  596 

116  958 

TYPE  PLATE  BENDIMG 
ELEHEBT  INCIDENCES 

201  421 

202  245 

203  1 5 2 

204  5 1 4 

205  532 

206  356 

207  263 

208  625 

209  754 

210  578 

211  485 

212  847 

213  865 

214  689 

215  596 

216  958 
TYPE  SPACE  FEAHE 
nENBEB  INCIDENCES 

5 1 10 

6 7 10 

7 2 11 

8 8 11 

9 3 12 

10  9 12 

11  10  11 
12  11  12 

13  MOA*  10 

14  12  •12A« 


ILEBEIT  GEIIRATIOV* 


TN-77-17(5-10) 


ELEMENT  PROPERTIES 

101  TO  116  TYPE  •CSTG'  THICK  .05 

201  TO  216  TYPE  'CPT'  THICK  .5313 

NENBER  PROPERTIES  PRISMATIC 

5 TO  12  AX  3.  IX  20.  lY  10.  IE  10. 

13  14  AX  30.  IX  200-  lY  100.  IE  100, 

CONSTANTS 

£ 10. E6  ALL 

POISSON  .3  ALL 

LOADING  1 

JOINT  DISPLACEMENTS 
• 10A*  DISP  X 1. 

LOADING  2 

JOINT  DISPLACEMENTS 
• 10A'  DISP  Y 1- 
LOADING  3 

JOINT  DISPLACEMENTS 
• 10A'  DISP  Z 1- 
LOADING  4 

JOINT  DISPLACEMENTS 
• 10A'  ROT  X 1 - 
LOADING  5 

JOINT  DISPLACEMENTS 
• 10A'  ROT  Y 1. 

LOADING  6 

JOINT  DISPLACEMENTS 
• 10A«  ROT  Z 1 . 

LOADING  7 

JOINT  DISPLACEMENTS 
•12A'  DISP  X 1. 

LOADING  8 

JOINT  DISPLACEMENTS 
•12A*  DISP  Y 1, 

LOADING  9 

JOINT  DISPLACEMENTS 
' 12A'  DISP  Z 1. 

LOADING  10 

JOINT  DISPLACEMENTS 

• 12A'  ROT  X 1 - 

LOADING  11 

JOINT  DISPLACEMENTS 

• 12A«  BOT  Y 1. 

LOADING  12 

JOINT  DISPLACEMENTS 

• 12A'  ROT  Z 1. 

DOHP  TIME 

STIFFNESS  ANALYSIS  REDOCE  BAND 

OOTPUT  BY  JOINT 

LIST  DISPLACEMENTS  ALL 

OOTPOT  DECIMAL  2 

OOTPOT  BY  MEMBER 

LIST  FORCES  13  14 

FINISH 


Fig. 5. 10-  Input  listing  of  STRUDL  jnodel  used  in  super  finite  element 
generation. 


53 


PNOSLBM 


SOPEB 


TITLE  - SKEWED  TBXAMGOLAE  80DEL 


TN-77-17(5-ll) 


NEHBEI  rOBCES 


HEHBEF 

13 


LOADING 

JOI  NT 

/ 

AXIAL 

SHEAR  T 

SHEAR  Z 

TORSION 

HOHENT  X 

HOHENT  Z 

1 

10A 

648215.69 

3.53 

105.98 

1578.36 

74398.31 

-478095.50 

10 

-648215.69 

-3.53 

-105.98 

-1578.36 

-74399.37 

478095.50 

2 

10A 

3.53'j 

L 14978.57 

802.67 

3185.35 

-21315.01 

374429.50 

10 

-3.53 

-14978.57 

-802.67 

-3185.35 

21306.98 

-374279.69 

3 

10A 

105.98  ^ 

L 802.67 

26598.18 

120934.06 

-665205.31 

18758.  24 

10 

-105.98  i 

f -802.67 

-26598.  18 

-120934.06 

664939.37 

-18750.21 

4 

10A 

1578.39 

1 3185.56 

1209  34.06 

4096715.00 

-3042881.00 

76  995,94 

10 

-1578,39  J 

' -3185.56 

-120934.06 

-4096715.00 

3041671.00 

-76964.06 

5 

10A 

74398. 31 

-21314.95 

-665205.31 

-3042880.00 

30660912.0 

-1263416.00 

10 

-74398.31 

21314.95 

665205.31 

3042880.00 

-30654256.0 

1263203.00 

6 

10A 

-4780  95.50 

374429.44 

18758. 19 

76994.81 

-1263416.00 

16246325.0 

10 

4780  95.50 

-374429.44 

-18758.  19 

-76994.81 

1263229.00 

-16242580.0 

7 

10A 

-648215.69 

-3.53 

-105.98 

-1578.39 

-74398.31 

478095.50 

10 

648215.69 

3.53 

105.98 

1578.39 

74399.37 

-478095.56 

8 

10A 

-3.53 

-14974.34 

-802.66 

-3185.53 

21314.95 

-374429.  12 

10 

3.53 

14974.34 

802.66 

3185.53 

-21306.92 

374279.37 

9 

10A 

-105.98 

-802.66 

-26581.77 

-120934.06 

665205.25 

-18758.15 

10 

105.98 

802.66 

26581.77 

120934,06 

-664939.37 

18750.12 

10 

1 OA 

-1578.39 

-3185.21 

-120934.06 

-4096714.00 

3042881.00 

-76993.19 

10 

1578.39 

3185.21 

120934.06 

4096714.00 

-3041672.00 

76961.37 

11 

10A 

-79699.25 

-18834.34 

-664415.06 

-3006241.00 

2612640.00 

325132.12 

10 

79699.25 

18834.  34 

664415.06 

3006241. OC 

-2605996.00 

^..^-325320.44 

2482628.00 

12 

10A 

47  82  72.00 

374587.37 

21391.12 

82340.31 

197239.44 

10 

-4782 72.00 

-374587.37 

-21391. 12 

-82340.31 

-197453.31 

-2478882.00 

PROBLEB  - 5DPE8 


TITLE  - STHBETBIC  TBIAIGOLAB  ROOEL 


HBH8EB  POBCES 


LOADING 

JOINT 

/ 

— //— 

HOHENTS 

AXXAL 

S8EAR  X 

SBEAB  Z 

TORSION 

HOHENT  X 

HOHENT  Z 

1 

10a 

648346. 87 

0.00 

0.00 

0.02 

78655,31 

-479487.50 

10 

-648346. 87 

-0.00 

-0.00 

-0.02 

-78655.31 

479487.50 

2 

10A 

0.00  > 

L 14986.55 

812.47 

3198.29 

-20319.92 

374886.81 

10 

-0.00  j 

1 -14986.55 

-812,47 

-3198,29 

20311.80 

-374736.94 

3 

10A 

0.00  1 

1 812.47 

26650.64 

121072.00 

-6  66327-25 

20319.88 

10 

-0.00  i 

f -812.47 

-26650.64 

-121072.00 

666060.69 

-20311.76 

4 

10A 

0.03 

\ 3198.21 

121072.06 

4110196.00 

-3028014.00 

79987.19 

10 

-0.03  J 

f -3198.21 

-121072.06 

-4110196.00 

3026803.00 

-79955.25 

5 

10A 

78655.  31 

-20319.94 

-666327.25 

-3028013.00 

30806912.0 

-1287286.00 

10 

-78655.31 

20319.94 

666327.25 

3028013.00 

-30800240.0 

1287083.00 

6 

10A 

-479487.50 

374886.75 

20319.84 

79986.69 

-1287286.00 

16272879.0 

10 

479487.50 

-374886.75 

-20319. 84 

-79986.69 

1287083.00 

-16269130.0 

7 

10A 

-648346.87 

-0.00 

-0.  00 

-0.02 

-78655.31 

479487.50 

10 

648346.87 

0.00 

0.00 

0.02 

78655.31 

-479487.50 

8 

10A 

0.00 

-14989.50 

-812.47 

-3198.18 

20319.95 

-374887.00 

10 

-0.00 

14989.50 

812.47 

3198.18 

-20311.83 

374737,12 

9 

10A 

0.00 

-812.47 

-26642.73 

-121074.81 

666327.69 

-20319,91 

10 

-0.00 

812.47 

26642.73 

121074.81 

-666061.25 

20311.79 

10 

10A 

0.02 

-3198.32 

-121074. 81 

-4110231.00 

3028016.00 

-79987.62 

10 

-0.02 

3198.32 

121074.81 

4110231.00 

-3026805.00 

79955.62 

11 

10A 

-78655.31 

-20319.86 

-666327.62 

-3028016.00 

2522800.00 

^270884.56 

10 

78655,31 

20319.86 

666327.62 

3028016.00 

-2516137.00 

27 108  7. 75 

12 

10A 

479487.56 

374887.00 

20319.95 

79987.37 

270884.69 

^ 2478964.00 

10 

-479487.56 

-374887.00 

-20319.95 

-79987.37 

-271087.87 

-2475215.00 

Figure  5.11. 


Comparative  STRUDL  results  for  differing  triangular  meshes 
(skewed  vs.  symmetric). 
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TYPE  SPACE  FRAME 


ELEMENT  INCIDENCES 
15  10  12 

ELEMENT  PROPERTIES 


TN-77-17(5-12) 


15  TYPE  'DUMMY'  THICK  1. 
STIFFNESS  MATRIX  GLOBAL 


DOF  DX  DY  DZ  RX  RY  RZ 
NODES  2 NDF  6 


SUBMATRIX  1 1 - 


648346. 

0. 

0. 

0. 

78655. 

-479487. 

0. 

14986. 

812. 

3198. 

-20319. 

374886. 

0. 

812. 

26650. 

121072. 

-666327. 

20319, 

0. 

3198. 

121072. 

4110196. 

-3028014. 

79987. 

78655. 

-20319. 

-666327. 

-3028013. 

30806912. 

-1287286. 

-479487. 

374886. 

20319. 

79986. 

-1287286. 

16272879. 

SOBH&THIX  2 
-648346. 

1 - 

-0. 

-0. 

-0. 

-78655. 

479487. 

0. 

-14989. 

-812. 

-3198. 

20319. 

-374887. 

0. 

-81  2. 

-26642. 

-121074. 

666327. 

-20319. 

0. 

-3198. 

-121074. 

-41 10231. 

3028016. 

-79987. 

-78655. 

-20319. 

-666327. 

-3028016. 

2522800. 

270884. 

479487. 

374887. 

20319. 

79987. 

270884. 

2478964. 

SUBMATRIX  2 
648346. 

2 - 

0. 

-0. 

-0. 

78655. 

-479487. 

0. 

14990. 

812. 

3197. 

20319. 

-374886. 

-0. 

81  2. 

26656. 

121072. 

666327. 

-20319. 

-0. 

3198. 

121072. 

4110195. 

3028013. 

-79985. 

78655. 

20319. 

666327. 

3028013. 

30806912. 

-1287286. 

-479487. 

-374886. 

-20319. 

-79985. 

- 1287286. 

16272876. 

Figure  5.12.  Partial  listing  of  STRUDL  super  finite  element  model 
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APPENDIX 


Description  of  Finite  Elements 
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Element  Type:  *BPR*  - Bending  Plate  Rectangle 


Shape  § Node  Order: 


Nodal  Unknowns:  DZ,  RX,  RY 


Results:  Planar  stress  couples 

and  strain  curvatures 
at  each  node. 


Element  Loads : 


Uniform  out-of-plane 
surface  forces 


References : 3 


Element  Characteristics:  BPR  is  a flat  rectangular  element  with  four  nodes, 
the  corner  points.  Its  transverse  displacement  expansion  is 

2 2 3 3 2 

w = + a^x  + a^y  + a^xy  + a^x  + a^y  + a^x  + a^y  + a^x  y 

2 3 3 

+ a^pxy  + a^^x  y + 


The  expansions  for  the  rotations  are  the  derivatives  of  that  for  w.  The  ele- 
ment satisfies  inter-element  compatibility  with  respect  to  translation,  but, 
since  only  two  nodal  values  for  normal  slope  are  available  along  each  side, 
the  quadratic  variation  is  not  unique  and  the  normal  slope  is  not  continuous. 
The  solution  will  converge  as  the  mesh  size  is  reduced,  but  the  convergence 
will  not  be  monotonic,  due  to  the  incompatibility  of  normal  slope. 
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Element  Type:  *CPT*  - Compatible  Plate  Triangle 


Shape  S Node  Order: 

Y, 


Nodal  Unknowns:  DZ,  RX,  RY 


Results;  Planar  stress  couples 
and  strain  curvatures 
at  each  node. 


Element  Loads: 


Uniform  out-of-plane 
surface  forces 


References : 


Element  Characteristics:  CPT  is  a flat  triangular  element  with  three  nodes, 
the  corner  points.  The  element  is  divided  into  two  triangles  with  the  fol- 
lowing displacement  expansion: 


2 2 

= ai  + a^x  + a^y  + a^x  + a^xy  + a^y 


3 2 3 

X + Ogxy  + a^y 


These  expansions  are  not  complete  cubics,  since  they  contain  no  x y terms. 

The  expansions  for  the  rotations  are  the  derivatives  of  that  for  w.  w and 
dw 

are  continuous  across  the  interior  boundary,  x = 0.  The  normal  slope 
varies  linearly  along  each  edge.  The  missing  terms  cause  displacement  com- 
patibility not  to  be  completely  satisfied  along  side  1-2,  since  the  displace- 
ment expression  involves  nodal  quantities  for  the  third  node.  The  accuracy 
of  the  element  is  quite  good,  however,  in  spite  of  this  incompatibility. 
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Element  Type:  *CSTG’  - Constant  Strain  Triangle 


Shape  § Node  Order: 

3 


Nodal  Unknowns:  DX,  DY 


Results:  Planar  stresses  and 

strains  at  centroid . 


2 


Element  Loads : 


In-plane  surface, 
body,  and  edge  forces 


References:  3,  5,  6 


Element  Characteristics:  CSTG  is  a triangular  element  with  three  nodes,  the 
corner  points.  The  element  stiffness  matrix  is  computed  based  on  displace- 
ment expansions  of  the  following  form 

X = + a^x  + a^y 

which  yield  a linear  displacement  field  over  the  element  and  along  each  side. 
The  stresses  and  strains  are  thus  constant  over  the  element. 
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Element  Type:  *IPCQ*  - Isoparametric  Cubic  Quadrilateral 


Nodal  Unknowns:  DX,  DY 

Results:  Planar  stresses  and 

strains  at  each  node. 


Element  Loads : 

In-plane  surface, 
body,  and  edge  forces 

References:  7,  8 


Element  Characteristics:  IPCQ  is  a flat  quadrilateral  element  with  12  nodes, 
the  4 corner  points,  and  the  8 third-points  of  each  side.  The  displacement 
expansions  are  of  the  form 

2 2 3 2 2 

X = + a^n  + + a^n  + + ci-jK  + OgC  n + OgCn 

3 3^  ,3 

+ a^pTi  + C + 

which  yields  a quartic  field  over  the  element,  but  a cubic  variation  along 
the  edges.  ^ and  n are  generalized  curvilinear  coordinates  which  vary  from 
-1  to  +1  across  the  element.  The  sides  of  IPCQ  may  be  initially  curved  if 
the  coordinates  of  the  midside  nodes  are  not  colinear  with  the  corner  nodes. 
If  coordinates  are  not  specified  for  the  midside  nodes,  they  will  be  com- 
puted assuming  that  the  sides  are  straight. 


Shape  S Node  Order: 


10  9 
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Element  Type:  *IPLQ*  - Isoparametric  Linear  Quadrilateral 

’IPLQCSH*  - Isoparametric  Linear  Quadrilateral  - Constant 

„ ■ Shear 

Nodal  Unknowns:  DX,  DY 


Results:  Planar  stresses  and 
strains  at  each  node. 


Element  Loads: 

In-plane  surface, 
body,  and  edge  forces 

References:  7,  8,  9 


Element  Characteristics:  IPLQ  is  a flat  quadrilateral  element  with  4 nodes, 
the  comer  points.  The  displacement  expansions  are  of  the  form 

X = + 0^5  + a^n  + 

which  yields  a quadratic  field  over  the  element,  but  a linear  variation  along 
the  edges.  ^ and  q are  generalized  coordinates  which  vary  from  -1  to  +1 
across  the  element. 

IPLQCSH  is  a flat  quadrilateral  element  with  4 nodes,  the  corner  points. 
The  displacement  expansions  are  the  same  as  the  IPLQ,  but  the  assumed  shear 
strain  variation  in  the  X-direction  is  different.  The  shear  strain  variation 
for  the  IPLQ  is  assumed  to  be  a constant  in  the  X-direction  and  equal  to  the 
centroidal  value  of  the  element;  whereas,  the  shear  strain  varies  linearly 
in  the  IPLQ.  The  resulting  formulation  tends  to  given  better  results  than 
the  IPLQ  for  problems  in  which  the  in-plane  forces  are  generated  by  overall 
bending. 

The  stiffness  matrix  is  generated  by  numerical  integration,  rather  than 
explicitly,  and  results  in  an  element  with  performance  characteristics  and 
accuracy  similar  to  the  PSRCSH  element.  The  formulations  are  identical  if 
the  element  is  rectangular  in  shape. 


Shape  § Node  Order: 

3 
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Element  Type:  ’IPLS’  - Isoparametric  Linear  Solid 

’IPLSCSH'  - Isoparametric  Linear  Solid  - Constant  Shear 


Nodal  Unknowns:  DX,  DY,  DZ 

Results:  Global  stresses  and 

strains  at  each  node. 


Element  Loads : 

Surface , body, 
and  edge  forces 


References:  7,  8,  9 

Element  Characteristics:  IPLS  is  a six-sided  solid  element  each  of  whose 
faces  is  quadrilateral.  The  element  has  8 nodes,  the  corner  points.  The 
displacement  expansions  are  of  the  form 

X = + a^n  + 

which  yields  a cubic  field  through  the  element,  a quadratic  field  on  the 
faces,  and  linear  displacements  along  the  edges.  ri>  ^ are  generalized 

coordinates  which  vary  from  -1  to  +1  across  the  element, 

IPLSCSH  is  a six-sided  solid  element  each  of  whose  faces  is  a quadrilat- 
eral. The  element  has  8 nodes,  the  corner  points.  The  displacement  expan- 
sions are  the  same  as  for  the  IPLS,  but  the  assumed  shear  strain  variation  in 
the  X-direction  is  different.  The  shear  strain  variation  for  the  IPLSCSU  is 
assumed  to  be  constant  in  the  X-direction  and  equal  to  the  centroidal  value 
of  the  element;  whereas,  the  shear  strain  varies  linearly  in  the  IPLS.  The 
resulting  formulation  tends  to  given  better  results  than  the  IPLS  for  the 
problems  in  which  the  structure  is  being  subjected  to  overall  bending. 


Shape  ^ Node  Order: 
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Element  Type:  *IPQQ*  - Isoparametric  Quadratic  Quadrilateral 


Nodal  Unknowns:  DX,  DY 

Results:  Planar  stresses  and 
strains  at  each  node. 


Element  Loads: 

In-plane  surface, 
body,  and  edge  forces 

References : 7,8 


Element  Characteristics:  IPQQ  is  a flat  quadrilateral  element  with  8 nodes, 
the  corner  points  and  the  midpoints  of  each  side.  The  displacement  expan- 
sions are  of  the  form 

2 2 2 2 

X = + a^n  + n + agCn 

which  yields  a cubic  field  over  the  element,  but  a quadratic  variation  along 
the  edges.  ^ and  n are  generalized  curvilinear  coordinates  which  vary  from 
-1  to  +1  across  the  element.  Although  the  stiffness  matrix  is  generated  by 
numerical  integration,  rather  than  explicitly,  the  result  is  an  element  with 
performance  characteristics  and  accuracy  similar  to  the  LSR  element.  The 
sides  of  IPQQ  may  be  initially  curved  if  the  coordinates  of  the  midside  nodes 
are  riot  colinear  with  the  corner  nodes.  If  coordinates  are  not  specified  for 
the  midside  nodes,  they  will  be  computed  assuming  that  the  sides  are  straight 


Shape  & Node  Order: 


7 3 
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Element  Type:  *IPQS*  - Isoparametric  Quadratic  Solid 


Shape  S Node  Order: 


Nodal  Unknowns:  DX,  DY,  DZ 


Results : 

Global  stresses 

and 

strains  at  each 

node. 

Element 

Loads : 

Surface,  body, 
and  edge  forces 

References:  7,  8 

Element  Characteristics:  IPQS  is  a six-sided  solid  element  each  of  whose 
faces  is  a quadrilateral.  The  element  has  20  nodes,  the  corner  points  and 
midpoints  of  each  edge.  The  displacement  expansions  are  of  the  form 

2 2 

X = + ctgCnC  + cCg^  + ot^^n 


2 2 2 

+ nc  + 

which  yields  a quartic  field  through  the  element,  a cubic  field  on  the  faces, 
and  quadratic  displacements  along  the  edges.  5,  q,  and  C are  generalized 
curvilinear  coordinates  which  vary  from  -1  to  +1  across  the  element.  The 
edges  of  IPQS  may  be  initially  curved  if  the  coordinates  of  the  midsidc  nodes 
are  not  colinear  with  the  corner  nodes.  If  coordinates  are  not  specified 
for  the  midside  nodes,  they  will  be  computed  assuming  that  the  edges  are 
straight . 
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Element  Type:  'LSR'  - Linear  Strain  Rectangle 


Shape  § Node  Order: 

u 


8 M 

' ‘ 1 

1 5 


Nodal  Unknowns:  DX,  DY 


3 


6 


Results:  Planar  stresses  and 

strains  at  the  corner 
nodes . 


Element  Loads: 


In-plane  surface, 
body,  and  edge  forces 


References : 10 


Element  Characteristics:  LSR  is  a flat  rectangular  element  with  8 nodes, 
the  corner  points  and  the  midpoints  of  each  side.  The  sides  must  be  straight 
initially.  The  displacement  expansions  are  of  the  form 

2 2 2 2 

X = + a2X  + a^y  + a^xy  + a^x  + a^y  + a^x  y + ^3^/ 

which  yields  a cubic  displacement  field  over  the  element,  but  a quadratic 
variation  along  the  boundaries.  The  displacement  along  each  boundary  depends 
only  upon  the  nodal  quantities  on  that  boundary,  satisfying  inter-element 
displacement  compatibility. 
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Element  Type:  *LST^  - Linear  Strain  Triangle 


Shape  5 Node  Order: 

3 


Nodal  Unknowns:  DX,  DY 


Results:  Planar  stresses  and 

strains  at  the 
comer  nodes. 


E lement  Loads ; 

In-plane  surface, 
body,  and  edge  forces 


References:  11 


Element  Characteristics:  LST  is  a flat  triangular  element  with  6 nodes,  the 
corner  points  and  the  midpoints  of  each  side.  The  sides  must  be  straight 
initially.  The  displacement  expansions  are  expressed  in  terms  of  the  tri- 
angular coordinates  C-j  follows: 


X = 


'(«i« 


2 ^ 


.i  = l 


These  expansions  yield  a complete  second-degree  polynomial,  giving  a quad- 
ratic displacement  field  over  the  element  and  along  each  boundary.  The  dis- 
placements along  each  boundary  depend  only  upon  the  nodal  quantities  on  that 
boundary,  satisfying  inter-element  displacement  compatibility. 
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Element  Type:  ’PBQl’  - Plate  Bending  Quadrilateral 


Shape  & Node  Order: 


4 


Nodal  Unknowns:  DZ,  RX,  RY 


Results : Planar  stress  couples 
and  strain  curvatures 
at  each  node. 


Element  Loads : 

Uniform  out-of-plane 
surface  forces 


References : 10 


Element  Characteristics:  PBQl  is  a flat  quadrilateral  element  with  4 nodes, 
the  corner  points.  The  stiffness  matrix  is  formed  as  the  average  of  the 
matrices  generated  by  summing  the  stiffnesses  of  the  two  pairs  of  triangles 
formed  by  diagonals  across  the  element. 


*^1234  ""  2 


The  triangles  used  in  this  case  are  CPT's.  The  resulting  solution  will  con- 
verge as  the  mesh  size  is  reduced,  and  the  results  are  more  accurate  than  the 
worst  combination  of  triangles,  but  less  accurate  than  the  best.  Being  com- 
posed of  CPT's,  the  absolute  accuracy  is  relatively  good  in  any  case.  The 
stress  results  are  best  when  the  triangular  areas  are  nearly  the  same;  i.e., 
the  quadrilateral  is  a parallelogram. 
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Element  Type:  ’PSR*  - Plane  Stress  Rectangle 

*PSRCSH*  - Plane  Stress  Rectangle  - Constant  Shear 


Shape  S Node  Order: 

Nodal  Unknowns:  DX,  DY 

Y. 

Results:  Planar  stresses  and 

strains  at  the  centroid. 

3 

4 

Element  Loads : 

b 

In-plane  surface. 

a 

body,  and  edge  forces 

1 

2 Xl 

References:  3,  5,  6,  9 

Element  Characteristics:  PSR  is  a rectangular  element  with  4 nodes,  the  cor- 
ner points.  The  element  stiffness  matrix  is  computed  based  upon  displacement 
expansions  in  the  local  coordinate  system  of  the  following  form 


X = + a^x  + a^y  + a^xy 

which  yields  a quadratic  displacement  field  over  the  element,  but  a linear 
variation  along  the  edges, 

PSRCSH  is  a rectangular  element  with  4 nodes,  the  corner  points.  The 
displacement  expansions  are  the  same  as  for  the  PSR,  but- the  assumption  con- 
cerning the  shear  strain  variation  in  the  X-direction  is  different.  The  shear 
strain  variation  for  the  PSRCSH  is  assumed  to  be  a constant  in  the  X-direction 
and  equal  to  the  centroidal  value  of  the  element;  whereas  the  shear  strain 
varies  linearly  in  the  PSR.  The  resulting  formulation  tends  to  give  better 
results  than  the  PSR  for  problems  in  which  the  in-plane  forces  are  generated 
by  overall  bending. 
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Element  Type:  *SBCT*  - Stretching  S Bending  Compatible  Triangle 


Nodal  Unknowns:  DX,  DY,  DZ,  RX,  RY 


Results:  Planar  stresses  and 
strains  at  centroid; 
planar  stress  couples 
and  strain  curvatures 
at  each  node. 


Element  Loads : 

Uniform  out-of-plane  surface 
forces.  In-plane  surface,  body, 
and  edge  forces . 


References:  12 

Element  Characteristics:  SBCT  is  a flat  triangular  element  with  3 nodes,  the 
corner  points.  The  element  stiffness  matrix  is  the  sum  of  the  CSTG  and  CPT 
stiffness  matrices  and  is  generated  as  an  overlay  of  those  two  elements.  The 
in-plane  and  out-of-plane  stiffnesses  and  displacements  are  uncoupled.  The  , 
in-plane  formulation  is  plane  stress. 


Shape  5 Node  Order: 
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Element  Type:  *TRIP*  - Triangular  Prism 


Shape  5 Node 

Order: 

Nodal  Unknowns:  DX,  DY,  DZ 

i 

Results:  Global  stresses  and 
strains  at  the  nodes. 

3 

/ ^ 

— ► V 

Element  Loads : 

1 A 

Surface,  body, 
and  edge  forces 

References : 

Element  Characteristics:  TRIP  is  a five-sided  solid  element  containing  two 


triangular  faces  and  three  quadrilateral  faces.  The  element  has  6 nodes,  the 
corner  points.  The  displacement  expansions  are  of  the  same  form  as  the  IPLS 
elements,  but  the  quadrilateral  front  and  back  faces  are  allowed  to  degenerate 
to  triangles.  This  expansion  yields  a quadratic  field  on  the  quadrilateral 
faces,  and  a linear  field  on  the  triangular  faces  and  along  the  edges. 
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Element  Type:  *UTLQ1*  - Unequal  Transition  Linear  Quadrilateral  - 1 Node 


Nodal  Unknowns:  DX,  DY 

Results:  Planar  stresses  and 

strains  at  the  nodes’. 


Element  Loads: 

In-plane  surface, 
body,  and  edge  forces 

References : 13 


Element  Characteristics:  UTLQl  is  a flat  triangular  element  with  4 nodes, 
the  corner  points  and  one  midpoint  on  the  first  side.  The  displacement  ex- 
plansions  are  expressed  in  terms  of  the  triangular  coordinates  ^2*  ^3 
as  follows: 

3 

X ' Z ('i  ■ ')  * ‘•^1^2 

i=l 

These  expansions  will  yield  a quadratic  displacement  field  over  the  element 
and  along  the  first  side,  and  a linear  displacement  field  over  the  remaining 
two  sides.  The  displacements  along  each  boundary  depend  only  upon  the  nodal 
quantities  on  that  boundary,  satisfying  inter-element  displacement  compati- 
bility. 


Shape  5 Node  Order: 

3 
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